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The effect of imperfections on surface critical properties is studied for Ising models with nearest- 
neighbour ferromagnetic couplings on simple cubic lattices. In particular, results of Monte Carlo 
simulations for flat, perfect surfaces are compared to those for flat surfaces with random, 'weak' 
or 'strong', interactions between neighbouring spins in the surface layer, and for surfaces with 
steps of monoatomic height. Surface critical exponents at the ordinary transition, in particular 
/3i = 0.80 ± 0.01, are found to be robust against these perturbations. 
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t^- ■ I. INTRODUCTION 

o\ : 

| Critical phenomena at perfect, flat surfaces of three-dimensional crystals have attracted much interest. In particular, 
. surface critical exponents have been estimated, both theoretically and experimentally. The general agreement is quite 
O ' satisfactory [1-3]. 

Nevertheless, there are still some open questions. For example, the critical exponent Pi of the surface magnetization 
at the ordinary transition in the Ising universality class seems to be known only to an accuracy of about ten percent, 
with f}\ ~ 0.80 ± 0.05, obtained from large scale Monte Carlo simulations [4,5], renormalization group calculations 
[6,7], and experiments on magnets and alloys such as FeCo [3,8]. 
' In addition, the role of surface imperfections, unavoidable in real materials, on critical properties has been studied 
. for three-dimensional systems in much less detail (for exceptions, see, e.g., mean-field theory [9] and renormalization 
C 1 group arguments [10] on surface layers with randomness), especially in simulations. 

G\ ', The aim of this article is twofold. Firstly, the values of critical surface exponents for perfect surfaces at the ordinary 
transition in Ising systems have been reexamined by computing effective exponents allowing to monitor easily the 
approach to the true asymptotic behaviour and to refine previous estimates on the asymptotic exponents, e.g., on the 
value of f3\. Secondly, we have also performed Monte Carlo simulations for two types of imperfections, corresponding 
to simple cases of amorphous and corrugated surfaces, see Fig. 1. The amorphous surface is mimicked by choosing 
randomly 'weak' or 'strong' nearest-neighbour ferromagnetic couplings between spins in the flat surface layer. A 
corrugation is introduced in a simple manner by placing a terrace of monoatomic height on half of the surface; the 
magnetization at the step edge is expected to reflect most strongly this form of perturbation as compared to the 
I 1 magnetization at the perfect surface. 
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II. ISING MODEL WITH A PERFECT SURFACE 

> ■ 

• i-H , 

To simulate the critical behaviour of the perfect surface, we study Ising models on simple cubic lattices with free 
boundary conditions for the top and bottom (surface) layers and periodic boundary conditions otherwise. The spins 
are denoted by S xyz (± 1), situated at lattice sites (x, y, z), where z numbers the layers, z = 1,2, M, with K x L 
spins per layer (x = 1,2, K; y — 1,2, L). Interactions are restricted to nearest-neighbours. Two different, 
ferromagnetic couplings may occur, depending on whether the neighbouring spins are at the surface [z = 1 or M), 
J s > 0, or not, Jb > 0, see Fig. la. 

The phase diagram of the corresponding semi-infinite (K,L,M — > oo ) Ising model is well established [1,2]. If the 
ratio of the surface coupling J s to the bulk coupling r = J s /Jb, is sufficiently weak, the system undergoes at 
the bulk critical temperature, ksT c / Jb ~ 4.5115 [11,12], an 'ordinary transition', with the bulk and surface ordering 
occurring at the same temperature, T c . Beyond a critical ratio, r > r c w 1.50 [5], the surface orders at a higher 
temperature, T s > T c , followed by the 'extraordinary transition' of the bulk at T c , see Fig. 2. At the critical ratio r c , 
one encounters the 'special transition point', with critical properties of the surface transition deviating from those at 
the ordinary or the distinct surface transition. 

The aim of our simulations has been to reexamine the surface critical properties at the ordinary transition, because 
of the rather wide spread of supposedly accurate and reliable estimates of the critical exponents, especially of the 
exponent of the surface magnetization, j3\. 

To study both surface and bulk properties, we computed various quantities, among others the profiles of the magne- 
tization per layer, m{z), and the susceptibility per layer, x( z )> defined by 
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For z — l(or M), one obtains the standard surface magnetization mi = m(z = 1) = m(z = M) and the susceptibility 
Xn = x( z = 1) = x( z = Af). xn describes the response of the surface magnetization to a surface field [1,2]. The 
response of mi to a bulk field is obtained from xi (which we also computed), given by [13] 
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(3) 



Most of the Monte Carlo simulations were performed using the one-cluster flip algorithm [14], augmented by runs 
using single spin flips. System sizes ranged from K x L x M = 25 x 25 x 20 to 150 x 150 x 120, thereby attempting 
to avoid finite-size effects on approach to the phase transition; i.e. computing properties of the semi-infinite system. 
To equilibrate the system and to extract thermal averages, we usually generated a few 10 5 clusters, and determined 
error bars from typically 5 realizations. Various ratios r = J s /Jb of the surface to the bulk couplings were chosen, in 
between r = and 2.0, studying predominantly the ordinary transition at r = 1.0. 

Typical magnetization profiles are depicted in Fig. 3, at r = 1.0 and T < T c . The magnetization increases from its 
surface value, mi, to the bulk value nib, at distances superceeding the bulk correlation length. To avoid finite-size 
effects, K = L and M have to be sufficiently large. Firstly, M may be chosen such that the profile m(z) displays a 
pronounced plateau around the center of the system. Monitoring then the if -dependence of mi, one may determine 
the suitable size of the Monte Carlo system. In Fig. 3, we included also the value of the bulk magnetization, as had 
been obtained in extensive and highly accurate previous simulations [12], thereby giving another check on the quality 
of the present data. 

In order to estimate the critical exponents of the magnetization, we define an 'effective exponent' (3 e ff(z,t) by 



Peff(z,t) = d\n(m(z))/d\n(t) 



(4) 



where t = \T — T c \/T c is the reduced temperature. Certainly, on approach to T c , t — ► 0, /3 e // becomes the asymptotic 
exponent f3(z), with (3{z = 1) = /3i being the critical exponent of the surface magnetization. 

Results for (3 e ff(z,t) are shown in Fig. 4, at fixed temperatures, t, and in Fig. 5, at fixed distances from the surface, 
z. Because the simulations are performed at discrete temperatures, U, (3 e ff is replaced by 



Peff(z,t) = ln(m(z,ti)/m(z,t i+1 ))/ln(ti/ti +1 ) 



(5) 



with t = (U + i;+i)/2. Very accurate Monte Carlo data are required to get reliable estimates. Error bars may 
be assigned to the effective exponents in various ways. We adopted the form 5/3 e ff = (\dm(z, ti)/m(z, U)\ + 
\5m(z,ti + i)/m(z,ti + i)\), where 5m(z,ti) follows from the variance of the magnetization data due to the different 
realizations, i.e. Monte Carlo runs with different random numbers. These error bars are rather conservative as com- 
pared, for instance, to those resulting from standard error propagation assumptions. 

From Fig. 5, one sees that the effective exponent of the surface magnetization (3 e ff(z = 1) increases almost linearly 
with t over a wide range of temperatures. Note that we omitted, for clarity, data which are obviously affected by 
finite-size effects. A linear extrapolation to the critical point yields fi\ = 0.80 ± 0.01. This value is expected to 
correct and refine the previous estimates [3-8], ranging from 0.75 [5] to 0.845 [7]. It is nicely consistent with the 
value estimated from a log-log plot for simulational data in between 0.01 < t < 0.1, with f}\ = 0.78 ± 0.02 [4] (being, 
presumably, slightly too low, because the rise of /3 e // with smaller t had been neglected). 

In the non-asymptotic region, the surface magnetization may be cast in the form mi w mo^' 3l (l + at x ). From the 
data depicted in Fig. 5, the effective exponent x of the corrections to scaling has been estimated to be close to one 
(x w 0.95), in the range 0.02 < t < 0.1. 
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Moving from the surface into the bulk, the effective exponent is usually lowered, see Fig. 4. Eventually, one may 
observe an interesting crossover phenomenon, as shown in Fig. 5. Away from criticality, T c , and sufficiently deep in 
the bulk (e.g. z = fO), the effective exponent follows closely the behaviour of the effective bulk critical exponent, 
as obtained from the corresponding magnetization data of the Ising model with full periodic boundary conditions 
[12]. However, on approach to T c , when the bulk correlation length becomes comparable to the distance from the 
surface (= 10 lattice spacings, in the example), /3 e // crosses over to its surface value. Indeed, at each distance (large, 
but finite) from the surface, one reaches, sufficiently close to T c , the surface exponent /3i « 0.80 instead of the bulk 
exponent (3 « 0.32. 

The susceptibility per layer, x( z )i displays close to T Cl at the ordinary transition, r — 1.0, a non-monotonic behaviour, 
see Fig. 6, with the maximal value shifting towards the center of the system as the temperature gets larger. In analogy 
to (5), one may define an effective exponent -f e ff(z). In the bulk, that exponent is clearly positive, with a value around 
one at t w 0.1, rising quite sharply to at least about 1.5 on further approach to T c (finally, finite-size effects tend 
to decrease the effective exponent again). However, \ remains, at T Cl finite at the surface, z = 1, with a cusp-like 
singularity, i.e. a negative critical exponent 711, by approaching criticality from high temperatures [1,2]. To estimate 
711, we determined the temperature T m at which, for a finite system (K,L,M), xn exhibits its maximum, xTi- A 
suitably defined effective exponent is 

(7ii)e//(*m) = -rfln(Xn - Xu{tm))/dln(t m ) (6) 

with t m = \T — T m \/T m . Obviously, (7n) e // — ► 711 as K,L,M — ► 00, T m — ► T c , and t m — > t — > 0. Based on 
our Monte Carlo data at r = 1.0, this ansatz leads to a rather rough estimate of 711, circumventing finite-size effects 
again, with 711 = —0.25 ± 0.1. This value may be compared to results of previous simulations [4], where only the sign 
of 711 was identified, and to predictions of renormalization group calculations [1,2,6,7], 711 « —0.33. 
The critical exponent for the response function of the surface magnetization to a bulk field, \\, may be estimated in 
a more standard and straightforward way, by computing in analogy to (5) the effective exponent 

(7i)e//(t) = -dln(Xi(*))/dln(t) (?) 

At r = 1, that exponent decreases with t, on approach to T c , approximately linearly in the range 0.02 < t < 0.1. A 
linear extrapolation yields the asymptotic exponent 71 = 0.78 ± 0.05, in agreement with previous estimates [1,2,4]. 
We performed additional, albeit much less extensive simulations near the special point, r = r c . In particular, at 
r = 1.50, (3 e ff(z = 1) is found to increase more strongly than linearly, at least down to t = 0.02, so that an 
extrapolation to t — ► is not obvious. A reasonable estimate seems to be f3\ = 0.23 ±0.01, in good agreement with a 
fairly recent large scale Monte Carlo study, fli — 0.237 ±0.005 [5], but a little bit lower than the field theoretical value, 
fii « 0.26 [7]. Similarly, our estimate for 71 (= 1.5 ± 0.1) is slightly higher than that suggested by a renormalization 
group calculation, 71 sa 1.30 [7]. In any event, the agreement is reasonably satisfactory. 

Beyond the special point, r > r c , the critical exponents at the surface transition, T s > T c , are believed to be in the 
universality class of the two-dimensional Ising model [1,2]. Indeed, our Monte Carlo data are consistent with that 
prediction.- At this point, we draw attention to recent Monte Carlo simulations on short-range correlation functions 
near the surface transition, T s , which have been performed to interpret spin-polarized photoelectron data of some 
magnets [15]. 

III. SURFACE WITH RANDOM COUPLINGS 

To mimic an amorphous surface, possibly due to (non-)magnetic impurities at the surface of the crystal, we replace 
the unique surface coupling J s by two random, ferromagnetic interactions, being either 'strong', J s i, or 'weak', J S 2, 
i.e. J s i > J S 2 > 0, see Fig. lb.- Note that an alternate form of randomness, namely random surface fields, had been 
considered before using Monte Carlo techniques [16]. 

Following detailed simulations of the two-dimensional dilute Ising model [17,18], both couplings were assumed to 
occur with the same probability. Then the ratio d = 3 s ij J s \ measures the degree of dilution: d = 1 corresponds 
to the perfect case, and d = to the percolation limit. We performed most of the simulations at d = 1/10, where 
critical dilution effects are expected to be easily detectable, because the crossover length to the dilution dominated 
critical regime is only a few lattice constants [17]. The ratio of the surface to bulk couplings r = (J s i + J S 2) I '(2J&) 
was varied in between 1.0 and 3.5, to study critical phenomena at the ordinary transition and to locate the special 
point, r = r c (d). For instance, r = 1.0 is realized when J si = (20/11) J;, and J s2 = (2/11) J;,, for d = 1/10. The sizes 
of Monte Carlo systems and lengths of runs were chosen like in the perfect case. 

Because the critical temperature (thence, the effective interaction) of the corresponding two-dimensional system is 
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reduced by randomness, at given mean coupling (J s i + J S 2)/2 [17], the special point may be argued to be shifted 
towards larger ratios r c , as dilution is increased (indeed, at d = 1/10, for example, we locate r c at 1.70 ± 0.1, see 
below, to be compared with 1.50 for the perfect surface). Accordingly, at r — 1.0 and d = 1/10 the ordinary transition 
is encountered. We shall discuss our findings for that case first, comparing it to results for the perfect surface, d = 1. 
Certainly, the transition temperature T c is not affected by the randomness in the surface couplings. 
At fixed temperature, T < T c , the dilution tends to decrease the magnetization m(z) at and near the surface, for 
distances small compared to the bulk correlation length. This trend is shown in Fig. 7 for the surface magnetization, 
mi. At first sight, perhaps, somewhat surprisingly, the energy per layer E{z) is also lowered compared to the perfect 
situation, but only in the surface layer. This behaviour may be explained by the ordering in clusters of strongly 
interacting, J s i(> J s ), spins, with a weak coupling between these, possibly oppositely oriented, clusters. 
Despite the obvious drop in m 1 due to the surface randomness, the ratio of m\{d = 1/10) /mi (d = 1) is almost 
constant, about 0.9, over a wide range of temperatures, from ksT ' j J& w 3.4 up to T c , with a very shallow minimum 
near ksT/Jb = 3.8. Therefore, the effective critical exponent e ff(z = 1), see Eq. (5), follows very closely that 
of the perfect case, as shown in Fig. 8, leading to the same estimate for f3\ — 0.80 ± 0.01, obtained from a linear 
extrapolation towards T c of the data for the effective exponent. Thence, the randomness in the surface couplings 
seems to be irrelevant for the asymptotic behaviour of the surface magnetization, and of minor importance even for 
the corrections to scaling. Similarly, the estimate for 71 is compatible with the one in the perfect case. Obviously, 
the numerical findings provide support to the conjecture by Diehl and Niisser, suggesting that short-range correlated 
dilution of surface interactions is irrelevant for surface critical exponents, based on a Harris-type criterion [19] and 
perturbation-theory to first-order [10].- The findings are corroborated by exact and simulational results on two- 
dimensional Ising models with surfaces, where (5\ is 1/2, independent of surface dilution [20,21].- Interestingly enough, 
in three-dimensional Ising models, random surface fields seem to be irrelevant as well [16]. 

To locate the special point, r — r c , one may determine the surface transition temperature, T s . At r < r c , T s — T c « 
4.5115 Jb/kB, while T s (r) > T c otherwise. Standard procedures may be applied, e.g., by analysing the finite-size 
dependence of the turning point in mi. The resulting phase diagram near the special point, for d = 1/10, is depicted 
in Fig. 9. 

However, close to the special point, T s deviates only minutely from T c , and ambiguities may arise. To get easily a lower 
bound for r c , one may monitor effective exponents, e.g., e ff(z = 1)- For instance, in the random case with d — 1/10, 
at r = 1.53, i.e. slightly above r c (d = 1), (3 e ff rises sharply at t < 0.05, being consistent with the asymptotic value 
/?i « 0.80. Even at r = 1.7, the effective exponent still tends to increase on approach to T c , to a value close to that 
of the perfect surface at the special transition, fli 0.23.- An estimate for an upper bound for r c may be obtained 
by assuming T s — T c , and looking for inconsistencies in the critical exponents. In particular, such an inconsistency 
occurs when /3 e ff(z = 1) becomes smaller than 1/8. Based on that kind of analysis, in addition to the standard 
procedures, we reached at an upper bound r c < 1.8. In summary, a reasonable estimate of r c (d — 1/10) seems to 
be 1.70 ± 0.1. To refine it, more sophisticated methods are needed [5], which, however, are beyond the scope of this 
study. The important point is, to notice that r c shifts towards a higher value when diluting the surface interactions. 



To study effects of corrugation on surface critical phenomena, we shall consider here thermal properties at a step 
of monoatomic height, superimposed on a perfect surface. 

As in the perfect case, we consider systems with K x L spins per (x, y) layer, with the exception of the topmost 
layer, z = 0, where spins are restricted to the center half, forming a strip-like terrace of monoatomic height along the 
y-direction, with two bordering steps, see Fig. lc. At the bottom, z = M, the surface is assumed to be flat. 
Three kinds of interactions may be introduced. Spins at the step-edges are coupled ferromagnetically, J e > 0. Other 
spins at the top (or bottom) of the system, having no upper (lower) neighbouring spins, interact through the surface 
coupling J s > 0, while all the remaining spins are coupled by Jf, > 0. The phase diagram in the (r — J s / Jb, ksT/ Jb)- 
plane is expected to be identical to that of the perfect surface [1], because the couplings at the one dimensional 
step-edges do not support additional long-range ordering. 

To quantify the influence of the additional terrace, we calculated the magnetization per row 



The magnetization at the step-edges, say, x — x s \,x S 2, mn = m(x s i,z = 0) = m(x s2 ,z = 0) is of particular 
interest, deviating possibly most significantly from the magnetization of the perfect surface, mi. Certainly, the aim 



IV. SURFACE WITH STEPS 
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is to identify the behaviour at a single step, i.e. the width of the terrace has to be large compared to the relevant 
correlation lengths, as it is the case in the limit of the semi-infinite system, (K,L,M — ► oo ). Therefore, careful 
finite -size analyses are needed again, considering the effect of all three linear dimensions, K, L, M, on the quantities 
of interest, especially m\\. For instance, on approach to criticality, T c , the length of the step, L, is of importance [21]. 
To monitor and, eventually, circumvent finite size effects, K ranged from 80 to 120, L from 80 to 240, and M from 
40 to 80. As before, we used the one-cluster-flip Monte Carlo algorithm, generating a few 10 5 clusters, and averaging 
typically over about 5 realizations. 

In Fig. 10, the magnetization per row, m(x, z), is depicted at fixed temperature near the ordinary transition, taking 
J e = J s = J;,, for various layers, starting at the top, z = 0, and moving into the bulk. The magnetization is minimal at 
the step-edge, due to the reduced coordination number at the step. To avoid finite-size effects, the system size has to 
be chosen in such a way, that, e.g., m(x, z) acquires (i) the numerically accurately known bulk magnetization [11,12] 
sufficiently far away from the surface, and (ii) the correct surface magnetization, mi, at the bottom, m(x,M), as 
determined in the perfect case. Furthermore, the influence of the length of the step, L, on the step edge magnetization 
needs to be scrutinized, thereby monitoring the impact of bulk and surface correlation lengths. 

The effective exponent of the step-edge magnetization, (/3 e ff)u(t), defined in complete analogy to (5), is shown in 
Fig. 1 1 . In the figure, Monte Carlo data affected by finite-size effects are included, characterised by a decrease in the 
effective exponent as one goes closer to T c . As observed readily, finite-size effects play an important role for unusually 
large systems even quite far away from T c . 

In the following we consider only the data unperturbed by finite-size effects. Comparing to Fig. 5, one notices an 
enhancement of the exponent (f3 e ff)n(t), at fixed distance from T c , i.e. at fixed t, relative to the effective exponent 
of the surface magnetization. In any event, (/3 e //)ii W seems to vary on approach to criticality almost linearly. Based 
on a linear extrapolation, the asymptotic value of 0n may be estimated to be 0.80 ± 0.015, in agreement with that 
of the magnetization of a flat surface. In contrast to the case of random couplings, corrections to scaling of the 
magnetization at the step-edge are distinctively different from those of mi for the perfect surface. 
Our result for j3n on surfaces with steps may be argued to be in accordance with the conjecture of Diehl and Niisser, 
stating that smooth corrugations correspond to irrelevant perturbations, in the sense of the renormalization group 
theory [10]. The conjecture, however, refers merely to the critical behaviour of the entire surface, and may not 
necessarily apply to the step-edge magnetization. 

V. SUMMARY 

Using Monte Carlo techniques, we studied critical properties at the ordinary transition, the special point, and 
the distinct surface transition of semi-infinite Ising systems on cubic lattices with nearest-neighbour ferromagnetic 
interactions. In particular, the influence of randomness and corrugations on the critical exponents at the ordinary 
transition has been investigated. 

Especially the critical exponent of the surface magnetization has been shown to be robust against both types of 
perturbations, as mimicked by random 'strong' and 'weak' couplings in the surface, and by steps of monoatomic 
height. Its value is always 0.80 ± 0.01(5). Even corrections to scaling are affected only rather mildly by the dilution, 
in contrast to the corrugated case comparing the step-edge magnetization to the surface magnetization, where the 
correction terms to the asymptotic power-law are found to differ appreciably. 

Our Monte Carlo findings support and extend conjectures of Diehl and Niisser, based on renormalization group ar- 
guments, suggesting the irrelevant character of these surface perturbations at the ordinary transition. The findings 
may also be helpful in interpreting experiments on surface critical phenomena, by showing that different 'dirt effects' 
are of minor importance for the asymptotic critical exponents. 
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Figure Captions 

Fig. 1: Geometry and interactions of the Ising models for (a) the perfect surface, and surfaces with (b) random couplings 
and (c) steps. 

Fig. 2: Schematic phase diagram for the semi-infinite three-dimensional Ising model; r denotes the ratio of the surface 
to the bulk couplings. 

Fig. 3: Magnetization profiles m(z) for the perfect surface with J s = Jb at, from bottom to top, fc^T/Jb = 4.49, 4.45, 
and 4.35, for systems of size 100 x 100 x 80. The dashed line denotes the bulk value [12]. 

Fig. 4: Effective exponent of the magnetization per layer, (3 e ff(z,t), of the Ising model with a perfect surface, with 
r = 1.0 at, from bottom to top, k B T/J b =4.25, 4.425, and 4.46, for systems of size 50 x 50 x 40. 

Fig. 5: (3 e ff(z,t) in the perfect case, with r = 1.0, at fixed distances from the surface, z, as a function of reduced 
temperature t, for systems of size 25 x 25 x 20 (triangles), 50 x 50 x 40 (diamonds), 100 x 100 x 80 (squares), 
and 150 x 150 x 120 (circles). The dashed line follows from bulk Monte Carlo data [12]. 
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Fig. 6: Susceptibility per layer \{ z ) m the perfect case with r — 1.0 at, from bottom to top, ksT/Jb = 4.2, 4.35, 4.45, 
and 4.47, for systems of size 50 x 50 x 40. 

Fig. 7: Surface magnetization mi vs. temperature, for the perfect and random surface, at r = 1.0, for systems of size 
50 x 50 x 40. 

Fig. 8: Effective exponent of the surface magnetization, (3 e ff{ z = 1), f° r a surface with random couplings, with J s2 = 
J sl /10 and r = 1.0. 

Fig. 9: Phase diagram for the three-dimensional Ising model with random couplings in the surface, J S 2 = J s \/10, close 
to the special point (at r c — 1.70 ± 0.1, as indicated in the figure). 

Fig. 10: Profiles of the magnetization per row m(x, z) for the Ising model with a terrace of monoatomic height on the 
surface, with J e = J s = J b , at T = 0.95T C , for systems of size 80 x 80 x 40. 

Fig. 11: Effective exponent of the step-edge magnetization, (/3 e //)n(i), as a function of reduced temperature t, with 
J e = J s = Jb, for various system sizes. 
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Pleimling/Selke: FIGURE 7 
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Pleimling/Selke: FIGURE 8 
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Pleimling/Selke: FIGURE 9 
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Pleimling/Selke: FIGURE 10 
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Pleimling/Selke: FIGURE 11 
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